Method of distributed estimation using multiple asynchronous sensors

ABSTRACT

A method and system are provided for merging data from a plurality of multiplexed measurement sources to a decision-maker. The method includes operations for receiving a corresponding plurality of measurements of the data, processing each measurement to respectively obtain local state estimates and local error covariances, determining a corresponding plurality of lag periods, offsetting each of the corresponding event times, supplying to a track fusion center the local state estimates and the local error covariances for summing the pluralities of the local state estimates as a fusion state estimate and the local error covariances as a fusion error covariance. The measurements to be fused are each acquired from its respective source and correspond to an associated sampling period within an acquisition interval. The lag periods represent a wait duration for obtaining the corresponding local state estimates and local error covariances. The system correspondingly includes a receiver to obtain the plurality of measurements, a set of associated processors to analyze the measurements and provide their local state estimates and the local error covariances; a time-offsetter to determine corresponding lag periods and displace each of the event times accordingly, and the track fusion center.

CROSS REFERENCE TO RELATED APPLICATION

Pursuant to 35 U.S.C. §119, the benefit of priority from provisional application 60/797,082, with a filing date of Apr. 28, 2006, is claimed for this non-provisional application.

STATEMENT OF GOVERNMENT INTEREST

The invention described was made in the performance of official duties by one or more employees of the Department of the Navy, and thus, the invention herein may be manufactured, used or licensed by or for the Government of the United States of America for governmental purposes without the payment of any royalties thereon or therefore.

BACKGROUND

Field of the Invention

The invention relates generally to a method for merging information about a target relative to a plurality of observers. In particular, the information represents relational characteristics such as position, velocity, acceleration based on range, azimuth and elevation with faster real-time response and reduced error compared to conventional techniques.

Threat assessment for decision-makers in the armed forces represents a vital element in detecting, assessing and engaging an approaching target. An observer may have one or more sensors (e.g., radar) that may have a different rate to determine entrance of an object within an observation perimeter. Upon such detection, a decision-maker may seek translational characteristics of the object to interpret the threat level posed therefrom. Such characteristics may include positional information of the object, such as range, azimuth and elevation relative to the sensor's vantage. Additional such characteristics may include change rate information such as speed and direction (i.e., velocity) of the approaching object at various levels of fidelity depending on the type of information sought. Identification of the object may include determining whether the object is friendly, neutral or hostile.

A conclusion of the object's hostility may necessitate a finer resolution (i.e., reduced uncertainty and time-delay) of positional relationships relative to the observers (i.e., sensors) available to the decision-maker. The more refined detail enables creation of a firing solution to engage the object, now declared a target. Such efforts to achieve such decisions and consequent engagement constitute tracking. The United States Navy employs the Aegis Combat System to coordinate data acquisition, interpretation and decision-making.

The Soviet Union's dissolution has led to a number of unforeseen consequences that are of concern to the United States military. The problem of particular concern to the United States Navy is the proliferation of inexpensive Theater Ballistic Missiles (TBMs) and Anti-Ship Ballistic Missiles (ASCMs) to unstable and hostile Third World countries. With troops likely to be deployed in such situations to defend both Allied and American interests, a rapidly evolving need exists for the Navy to operate effectively in littoral environment. The Navy has been designed to operate in a blue-water environment against a Soviet-style threat. Naval adaptability is such that there can be a redeployment to the new environment provided that there is some sharing of data among the multiple sensor/weapons platforms (e.g., Aegis-equipped platforms, such as cruisers and destroyers primarily).

SUMMARY

Tracking with multiple sensors has long been recognized as an effective method to obtain an accurate picture of a surveillance area of interest. The use of multiple sensors poses challenging problems in term of communications and processing. To generate a common picture at the different nodes, all the sensor data must be processed simultaneously by all of the nodes. Performing sequential processing represents the optimal method of processing the data. However, this requires that all the nodes receive the same data, at the same time, and without delays. A real-time optimal solution to the tracking problem may not be feasible due to the large amount of data to be distributed, the different variety of sensors, and the limited bandwidth of current communications systems.

Conventional multiple-sensor tracking yield disadvantages addressed by various exemplary embodiments of the present invention. In particular, various exemplary embodiments provide a method for merging data from a plurality of multiplexed measurement sources to a decision-maker. The method includes operations for receiving a corresponding plurality of measurements of the data, processing each measurement to respectively obtain local state estimates and local error covariances, determining a corresponding plurality of lag periods, offsetting each of the corresponding event times, supplying to a track fusion center the local state estimates and the local error covariances for summing the pluralities of the local state estimates as a fusion state estimate and the local error covariances as a fusion error covariance. The measurements to be fused are each acquired from its respective source and correspond to an associated sampling period within an acquisition interval. The lag periods represent a wait duration for obtaining the corresponding local state estimates and local error covariances.

Various exemplary embodiments provide a system that correspondingly includes a receiver to obtain the plurality of measurements, a set of associated processors to analyze the measurements and provide their local state estimates and the local error covariances, a time-offsetter to determine corresponding lag periods and displace each of the event times accordingly, and the track fusion center.

In various exemplary embodiments include calculating a corresponding plurality of weighting factors and multiplying them to each respective local state estimate, and then summing the weighted state estimates for the fusion state estimate. Various additional embodiments include summing the local error covariances as the fusion error covariance by minimizing a mean square error term.

BRIEF DESCRIPTION OF THE DRAWINGS

These and various other features and aspects of various exemplary embodiments will be readily understood with reference to the following detailed description taken in conjunction with the accompanying drawings, in which like or similar numbers are used throughout, and in which:

FIG. 1 is a block diagram view of an asynchronous track fusion architecture.

DETAILED DESCRIPTION

In the following detailed description of exemplary embodiments of the invention, reference is made to the accompanying drawings that form a part hereof, and in which is shown by way of illustration specific exemplary embodiments in which the invention may be practiced. These embodiments are described in sufficient detail to enable those skilled in the art to practice the invention. Other embodiments may be utilized, and logical, mechanical, and other changes may be made without departing from the spirit or scope of the present invention. The following detailed description is, therefore, not to be taken in a limiting sense, and the scope of the present invention is defined only by the appended claims.

This disclosure is organized in the order of: (i) introduction; (ii) advantages of the inventive approach; (iii) description of the asynchronous track fusion problem; (iv) calculation of weighted matrices to express the covariance matrix; (v) an asynchronous multisensor track fusion filter that yields an analytical minimum mean square solution for an arbitrary number of asynchronous tracks; (vi) a description of the architectural diagram with timeline to provide the solution; (vii) application for tracking multiple targets with multiple sensors; (viii) derivation of the Bar-Shalom—Campo fusion rule as a simplified example case of the fusion algorithm; (ix) concluding remarks; and (x) appendices for proofs.

Introduction: New processing techniques need to be pursued in order to effectively take advantage of the abundance of multiple sensor data while considering such issues as communication delays, data arriving out of sequence, the bursty nature of the data (i.e., sudden arrival of information bursts before and after long stasis intervals). Such filtering techniques are known as asynchronous multi-sensor fusion approaches. With data sharing among the sensor platforms, many of the problems associated with the littoral environment can be mitigated by well-designed algorithms. The most important use of such data enables an effective picture of the battle-space environment to be formed that minimizes the problems of miscorrelations between threats and friendlies, creation of false tracks, increased battle-space awareness, and improved reaction time against high speed TBMs and low observable targets, such as ASCMs.

The idealized single integrated air picture (SIAP) is described as a common picture of the threat environment from the viewpoint of an omniscient observer (God's viewpoint) shared to the widest possible extent (or even cloned if possible) among every sensor and weapons node that constitutes the SIAP network. This common tactical picture would be displayed to decision-makers as well as shooters at all levels across the node. This becomes a difficult problem to solve for real-time applications when there is a large amount of data being distributed across the network. Additional difficulties occur when the sensors are dissimilar and with the limited bandwidth of current communication systems.

For an idealized “system-of-the-combat-system” the data from the individual sensors may be propagated to processing node(s) on a variety tactical platforms by a communication network. For each node a, there is a node b, such that b=1, 2, . . . , N that may be connected to node a. A message sent from node a to a specific element, say M, arrives at a time which differs from its arrival time at N, so that t_(aN)≠t_(aN). Furthermore, the network being temporally associative is not always a valid assumption, so t_(aN)≠t_(aM). Because of this problem, the data are not well-ordered. Merely because the data from sensor a are received before the data from sensor b, does not imply that the measurement or track data from sensor a are “earlier” than sensor data from sensor b. The data are not inherently time ordered. Because the sensors can have different data rates, they measure the target position and/or velocity at different times. Hence the data are asynchronous. The rate of data arrival not being exactly periodic represents another adverse temporal feature of the data. Real-time systems are only approximately periodic, so the variation in the rate of arrival resembles a Poisson process about the mean. (Some literature refers to this as “bursty” message traffic, i.e., subject to episodic eruption of events in an otherwise stasis background.) The temporal nature of the message traffic across shared data networks may be best described in general as delayed, Poisson, and asynchronous.

Most attempts to generate a composite track have been based on the assumption that the data are ordered sequentially and that an estimation algorithm such as the Kalman filter or the interactive multiple model (IMM) filter can be applied to process the data as they arrive at the processing node in real time. Given the actual nature of the data, this poses a significant difficulty because in most cases the data that arrive from individual sensors are not temporally ordered, which eliminates the possibility of optimally processing the data. Asynchronous data can produce a tactical picture that is inconsistent with reality when compared to properly ordered data. An alternative to the batch processing approach is to window the data within a fixed time interval and process all the data in a manner that allows for optimal processing. That process involves a first step to properly time order the data.

Conventional fusion methods use a single processor to process all of the different sensor data sequentially from a given target. The proposed multi-sensor fusion approach uses a two-level processing scheme comprised of local and global processing levels. At the local level, each processing node processes its local data in real-time and produces a local track at the time when the local data were acquired. These local tracks are then communicated, at different times, to a fusion center for further processing. The fusion process waits for an appropriate and/or adaptive time before fusing all of the incoming tracks into a global track for a given target. The time between updates at the fusion center depends on such factors as the computational loading of the processors, the communication delays between local processing nodes and the track fusion node, the type of motion of the target, etc.

The use of multiple sensors for state estimation can lead to better quality estimates than when a single sensor is used. For the same number of resources, sensors with different data rates (i.e., asynchronous sensors) can provide better coverage than synchronous sensors. The local estimates produced by asynchronous sensors are not time coincident. In the, presence of communication delays between sensor platforms and fusion center and in the presence of out-of-sequence data, the optimal centralized processing architecture may not be feasible in real time.

Taking into consideration the computational load sharing and the survivability of the system, the distributed estimation approach becomes an attractive alternative. On the other hand, current distributed estimation architectures assume that the sensors used are synchronous and no communication delays exist between the fusion center and sensor platforms. Many of the existing track fusion algorithms are optimal only for synchronized tracks. In order to perform synchronous track fusion when practical asynchronous sensors are used, the local estimates have to be synchronized at a considerable additional processing and computational cost.

Track fusion (also known as distributed filtering) represents the merging of state-estimates or tracks from multiple sensors, distinguishable either by type of measurement acquired or vantage relative to target. The information acquired thereby can vary according to several parameters, such as signal-to-noise quality, range resolution, azimuth and elevation resolutions, distinction over clutter, etc.

However, out-of-sequence tracks have not been considered in the literature when dealing with distributed estimation. The advantages from various embodiments of this invention include solving a general sensor-to-sensor track fusion problem that is directly applicable to an arbitrary number of asynchronous tracks where communication delays exist between sensor platforms and fusion center and where some of the tracks arrive out-of-sequence. The track fusion solution is optimal in the minimum mean square sense for the proposed linear fusion rule. No special processing, such as prediction or smoothing is needed in the proposed solution to synchronize the local tracks. In addition, using the proposed batch processing approach, latent and out-of-sequence tracks can be included without additional processing and computational cost.

Advantages: The benefit to combining data generated from single sensor platforms to synthesize a common picture generated from the data from all the sensor platforms to create a more complete picture of the underlying operational environment are legion. A more complete and accurate picture that one can base actionable decisions upon, e.g., those decisions that lead to identification and characterization of a tracked object as well as deployment of weapons or other resources with high probability of success, has many advantages to the war-fighter. These advantages include:

-   -   Reaction time is improved because the threat is detected         earlier.     -   A more unified understanding of the threat environment is         available, so resources to engage threats can be better         allocated.     -   Better state estimates are available on threats, so their         behavior can be better predicted, and therefore optimal         engagement strategies developed.     -   Environmental penetration is improved, due to different look         angles at the same thereat, so favorable geometry can be taken         advantage of by the war-fighter.     -   The ability to lessen point failures due to reliance on single         platform rather than multi-platform performance.

These are only some of the benefits of a unified operational environment that can only be synthesized from single platforms that share data to produce a common picture. This picture can be produced from either a track fusion or a data fusion approach to data gathered from single sensor platforms. The data fusion approach is too computational intensive and requires too much communication bandwidth for near term deployment on sensor platforms, so the track fusion approach is preferred. The proposed asynchronous track fusion deals with all the practical problems associated with track fusion in a real world environment that treats all of the difficulties with real world track fusion.

The Asynchronous Multi-Sensor Track Fusion Problem: Consider a dynamical system, such as a moving target, whose state evolves according to the following linear stochastic differential equation {dot over (X)}(t)=A X(t)+G W (t),  (2.1) where X(t) is the continuous-time system state at time t and Ŵ(t) is a zero mean white Gaussian process noise with covariance defined as E[ W (t) W(τ)′]= q (t)δ(t−τ),  (2.2) and the prime symbol ′ denotes transpose.

Assume that the state of the system is observed by a finite but arbitrary number, m, of sensors that have different data rates. Let the discrete time measurement z_(i) taken by sensor i, (i=1, 2, . . . , m) at time t_(ki) (within an interval k), be given by z _(i)(t _(k) _(i) )=H _(i)(t _(k) _(i) )X(t _(k) _(i) )+V _(i)(t _(k) _(i) ),  (2.3) where

$\begin{matrix} {{X\left( t_{k_{i}} \right)} = {{{\Phi\left( {t_{k_{i}},t_{k_{i} - 1}} \right)}{X\left( t_{k_{i} - 1} \right)}} + W_{t_{k_{i} - 1}}^{t_{k_{i}}}}} & (2.4) \\ {{\Phi\left( {t_{k_{i}},t_{k_{i} - 1}} \right)} = e^{A{({t_{k_{i}} - t_{k_{i} - 1}})}}} & (2.5) \\ {W_{t_{k_{i} - 1}}^{t_{k_{i}}} = {\int_{k_{i} - 1}^{k_{i}}{{\Phi\left( \ {t_{k_{i}},\tau} \right)}G{\overset{\_}{W}(\tau)}{\mathbb{d}\tau}}}} & (2.6) \\ {{E\left\lbrack {W_{t_{k_{i} - 1}}^{t_{k_{i}}}\left( W_{t_{k_{i} - 1}}^{t_{k_{i}}} \right)}^{\prime} \right\rbrack} = Q_{t_{k_{i} - 1}}^{t_{k_{i}}}} & (2.7) \\ {{Q_{t_{k_{i} - 1}}^{t_{k_{i}}} = {\int_{k_{i} - 1}^{k_{i}}{{\Phi\left( \ {t_{k_{i}},\tau} \right)}G{\overset{.}{q}(\tau)}G^{\prime}{\Phi\left( {t_{k_{i}},\tau} \right)}^{\prime}{\mathbb{d}\tau}}}},} & (2.8) \end{matrix}$ such that X(t_(ki)) is the system state at time t_(ki), φ(t_(ki), t_(ki-1)) is the transition matrix, W represents the state process noise, E represents the expected value, and Q represents the process noise covariance matrix.

The measurement noise V_(i) (t_(ki)) is a zero mean white Gaussian process with

$\begin{matrix} {{{E\left\lbrack {{V_{i}\left( t_{k_{i}} \right)}{V_{i}\left( t_{k_{m}} \right)}^{\prime}} \right\rbrack} = {{R_{i}\left( t_{k_{i}} \right)}\delta_{t_{k_{i}},t_{k_{i}}}}},} & (2.9) \end{matrix}$ where t_(ki)=k_(i)T_(i), and T_(i) is the sampling period of sensor i. Moreover, R_(i)(t_(ki)) represents sensor i noise covariance and δ is the Kronecker delta with respect to the time interval.

In target tracking applications, the system dynamics is linear in the Cartesian coordinates while the measurements are nonlinear function of the state. In such cases, the measurement model may be obtained by linearization around the state estimate.

Each sensor platform may be equipped with computational capability to maintain its own local track by processing the sensor measurements from eqn. (2.3) as they are received using the system dynamical model from eqn. (2.4). The local state estimate produced by processing the measurement z_(i) (t_(ki)) of sensor i, taken at time t_(ki), is denoted by X_(i)(t_(ki)|t_(ki)). The Kalman filter update of such an estimate is given by

$\begin{matrix} {{{X_{i}\left( {t_{k_{i}}❘t_{k_{i}}} \right)} = {{{A_{i}\left( t_{k_{i}} \right)}{\Phi\left( {t_{k_{i}},t_{k_{i} - 1}} \right)}{X\left( {t_{k_{i} - 1}❘t_{k_{i} - 1}} \right)}} + {{K_{i}\left( t_{k_{i}} \right)}{V_{i}\left( t_{k_{i}} \right)}} + {{K_{i}\left( t_{k_{i}} \right)}{{H_{i}\left( t_{k_{i}} \right)}\left\lbrack {{{\Phi\left( {t_{k_{i}},t_{k_{i} - 1}} \right)}{X_{i}\left( t_{k_{i}} \right)}} + W_{t_{k_{i} - 1}}^{t_{k_{i}}}} \right\rbrack}}}},} & (2.10) \end{matrix}$ where A_(i)(t_(k) _(i) )=[I−K_(i)(t_(k) _(i) )H_(i)(t_(k) _(i) )], such that I is the identity matrix, K_(i)(t_(ki)) is the gain of the filter used to process the measurements acquired by platform i and H_(i)(t_(ki)) is the measurement matrix of sensor i.

The local state estimate error associated with sensor i at time t_(ki) can be defined as a difference between the local state estimate and the system state. {tilde over (X)} _(i)(t _(k) _(i) |t _(k) _(i) )=X _(i)(t _(k) _(i) |t _(k) _(i) )−X _(i)(t _(k) _(i) ).  (2.11)

Substituting eqns. (2.4) and (2.6) into (2.11) enables the error in the estimate of local track i to be expressed as

$\begin{matrix} {{{\overset{\sim}{X}}_{i}\left( {t_{k_{i}}❘t_{k_{i}}} \right)} = {{{A\left( t_{k_{i}} \right)}\left\lbrack {{{\Phi\left( {t_{k_{i}},t_{k_{i} - 1}} \right)}{\overset{\sim}{X}\left( {t_{k_{i} - 1}❘t_{k_{i} - 1}} \right)}} - W_{t_{k_{i} - 1}}^{t_{k_{i}}}} \right\rbrack} + {{K\left( t_{k_{i}} \right)}{{V\left( t_{k_{i}} \right)}.}}}} & (2.12) \end{matrix}$

Depending on the individual data rate of the sensors within the fusion network, a sensor may obtain several measurements during a given time interval, while another sensor obtains none or only one measurement. The associated local processor updates the local track as many times as the number of measurements taken during that interval.

Calculation of Weighted Matrices: In order to improve the quality of the state estimate, local tracks are communicated from the sensors' platform to a central site for the purpose of estimation fusion. Various exemplary embodiments apply to an arbitrary communication rate. The sensor platforms can communicate their track after every update, after a given number of updates, or after every time period. Communication delays may exist between the sensor platforms and the fusion center. As an example, t=t_(k-1), represents the last time the fusion center performed a track fusion, while t_(k)=t_(k-1)+T_(f) represents the time of the next fusion time.

The period T_(f) represents an adaptive design parameter. This period can be as small as the time necessary to receive two tracks from two different sensor platforms or as long as the time necessary to receive as many tracks as the total number of sensors in the network. The value of T_(f) depends primarily on the priority given to the track fusion task, which depend on the application being processed. For example, in target tracking applications when a hostile target presents a high level of threat, the waiting time T_(f) should be small compared to the case where the target is in a constant velocity mode.

The fusion rule in the described exemplary embodiments is applicable regardless of the value of T_(f). For a case such that during the time interval [t_(k-1), t_(k)], a number of asynchronous tracks arrive at the fusion center. Because of communication delays, some of the tracks are generated during the time interval [t_(k-1), t_(k)] and some may have been generated during the interval [t_(k-2), t_(k-1)] and arrive during the interval [t_(k-1), t_(k)]. Before the track fusion process takes place, some of the tracks may be removed. If more than one track generated by the same sensor is received, only the most recent track is used for fusion.

If only one track is received from a given sensor platform and was obtained during the time interval [t_(k-2), t_(k-1)], a decision has to be made whether to keep such track or not. This mainly depends on the data rate of the sensor, the sensor accuracy, as well as the respective age of that track. The tracks that will be kept for track fusion are called valid tracks.

Therefore, the number n of tracks to be fused during the interval [t_(k-i), t_(k)] is an arbitrary number between 2 and m, where m is the number of sensors used. Because of the difference in communication delays, some of the validated tracks may arrive out-of-sequence. To deal with the out-of-sequence tracks without special processing, all the tracks that arrive at the fusion center during the time [t_(k-1), t_(k)] after removing redundant tracks, are fused simultaneously. The asynchronous multi-sensor track fusion problem can be stated as follows:

Given a number of asynchronous valid tracks, {X_(i)(t_(k)|t_(k)), P_(i)(t_(k)|t_(k))}, i=1, 2, . . . , n, that arrive at the fusion center during the time interval [t_(k-i), t_(k)], find the best estimate in the minimum mean square sense of the system state at time t_(k) when computed according to the fusion rule in eqn. (2.13).

$\begin{matrix} {{{X_{f}\left( t_{k_{i}} \middle| t_{k_{i}} \right)} = {\sum\limits_{i = 1}^{n}\;{{L_{i}\left( t_{k_{i}} \right)}{X_{i}\left( t_{k_{i}} \middle| t_{k_{i}} \right)}}}},} & (2.13) \end{matrix}$ where the L_(i)'s are unknown matrices to be determined.

The error of the fused track at time t_(k)=k T_(f) can be defined as {tilde over (X)} _(f)(t_(k)|t_(k))=X _(f)(t_(k)|t_(k))—X(t _(k)).  (2.14)

The equations for system dynamics (2.4) and the fusion rule (2.13) yields

$\begin{matrix} {{{\overset{\sim}{X}}_{f}\left( t_{k} \middle| t_{k} \right)} = {{\sum\limits_{i = 1}^{n}\;{L_{i}{X_{i}\left( t_{k_{i}} \middle| t_{k_{i}} \right)}}} + {\left\lbrack {{\sum\limits_{i = 1}^{n}\;{L_{i}{\Phi\left( {t_{k_{i}},t_{k}} \right)}}} - I} \right\rbrack{X\left( t_{k} \right)}} + {- {\sum\limits_{i = 1}^{n}\;{L_{i}{\Phi\left( {t_{k_{i}},t_{k}} \right)}{W_{t_{k_{i}}}^{t_{k}}.}}}}}} & (2.15) \end{matrix}$

If all the local filters are unbiased, then second term of (2.15) is zero and the fused track is unbiased if the sum over all the tracks of weighted estimates is equal to the identity matrix:

$\begin{matrix} {{\sum\limits_{i = 1}^{n}\;{L_{i}{\Phi\left( {t_{k_{i}},t_{k}} \right)}}} = {I.}} & (2.16) \end{matrix}$

Eqn. (2.16) represents the first constraint on the choice of the weighting matrices. Solving for L_(n) by using (2.16) yields the second constraint equation:

$\begin{matrix} {L_{n} = {{\Phi\left( {t_{k_{i}},t_{k_{n}}} \right)} - {\sum\limits_{i = 1}^{n}\;{L_{i}{{\Phi\left( {t_{k_{i}},t_{k_{n}}} \right)}.}}}}} & (2.17) \end{matrix}$

Inserting eqns. (2.16) and (2.17) into (2.15) allows rewriting the error in the fused track at time t_(k) as

$\begin{matrix} {{{\overset{\sim}{X}}_{f}\left( t_{k} \middle| t_{k} \right)} = {{\sum\limits_{i = 1}^{n}\;{L\left\lbrack {{\overset{\sim}{X}\left( t_{k_{i}} \middle| t_{k_{i}} \right)} - {{\Phi\left( {t_{k_{i}},t_{k_{n}}} \right)}{\overset{\sim}{X}\left( t_{k_{n}} \middle| t_{k_{n}} \right)}} - {{\Phi\left( {t_{k_{i}},t_{k}} \right)}\left\{ {W_{t_{k_{i}}}^{t_{k}} - W_{t_{k_{n}}}^{t_{k}}} \right\}}} \right\rbrack}} + {{\Phi\left( {t_{k_{i}},t_{k_{n}}} \right)}{\overset{\sim}{X}\left( t_{k_{n}} \middle| t_{k_{n}} \right)}} - {W_{t_{k_{n}}}^{t_{k}}.}}} & (2.18) \end{matrix}$

The error covariance matrix of the fused track P_(f) can consequently be defined as

$\begin{matrix} {{{P_{f}\left( {t_{k}❘t_{k}} \right)} = {E\left\{ {{{\overset{\sim}{X}}_{g}\left( {t_{k}❘t_{k}} \right)}{{\overset{\sim}{X}}_{f}\left( {t_{k}❘t_{k}} \right)}^{\prime}} \right\}}},} & (2.19) \end{matrix}$ enabling determination of the weights L_(i) that define the optimal asynchronous track fusion filter.

The Asynchronous Multisensor Track Fusion Filter: There are two steps to solving the asynchronous track fusion problem outlined in the previous section. First, the covariance matrix of the fused track may be expressed as a function of the weight matrices L_(i) using the fusion rule in eqn. (2.13). Second, the matrices may be analytically determined by minimizing the mean square error of the fused track.

For 1≦i<n−1 and 1≦j≦n−1, define M _(ij) =P _(ij)(t _(k) _(i) ,t _(k) _(j) )+φ(t _(k) _(i) ,t _(k) _(n) ) [P _(nn)(t _(k) _(n) ,t _(k) _(n) )+Q _(w)(t _(k) _(i) ,t _(k) _(j) ,t _(k) _(n) )]φ(t _(k) _(j) ,t _(k) _(n) )′−P _(in)(t _(k) _(i) ,t _(k) _(n))φ(t _(k) _(j) ,t _(k) _(n) )′−φ(t _(k) _(i) ,t _(k) _(n) )P _(nj)(t _(k) _(n) ,t _(k) _(j) )−P _(xw)(t _(k) _(i) ,t _(k) _(j) )φ(t _(k) _(j) ,t _(k))′−φ(t _(k) _(i) ,t _(k))P _(xw)(t _(k) _(j) ,t _(k) _(i) )′+φ(t _(k) _(i) ,t _(k) _(n) )P _(xw)(t _(k) _(n) ,t _(k) _(j) )φ(t _(k) _(j) ,t _(k))′+φ(t _(k) _(i) ,t _(k))P _(xw)(t _(k) _(n) ,t _(k) _(i) )′φ(t _(k) _(j) ,t _(k) _(n) )′  (3.1)

$\begin{matrix} {M_{n} = {{{\Phi\left( {t_{k},t_{k_{n}}} \right)}{P\left( {t_{k_{n}},t_{k_{n}}} \right)}{\Phi\left( {t_{k},t_{k_{n}}} \right)}^{\prime}} + Q_{t_{k_{n}}}^{t_{k}}}} & (3.2) \\ {N_{i} = {{{P_{in}\left( {t_{k_{i}},t_{k_{n}}} \right)}{\Phi\left( {t_{k},t_{k_{n}}} \right)}^{\prime}} - {{{\Phi\left( {t_{k_{i}},t_{k_{n}}} \right)}\left\lbrack {{P_{nn}\left( {t_{k_{n}},t_{k_{n}}} \right)} - Q_{t_{k_{i}}}^{t_{k_{n}}}} \right\rbrack}{\Phi\left( {t_{k},t_{k_{n}}} \right)}^{\prime}}}} & (3.3) \end{matrix}$ P_(ij)(t _(k) _(i) ,t _(k) _(j) )=A(t)[φ(t _(k) _(i) ,t _(k) _(i) ₋₁)P(t _(k) _(i) ₋₁ ,t _(k) _(j) ₋₁)φ(t _(k) _(j) ,t _(k) _(j) ₋₁)′+Q _(c)(t _(k) _(i) ,t _(k) _(j) )]A _(j)(t _(k) _(j) )′+K _(i)(t _(k) _(i) )R _(i)(t _(k) _(i) )K _(j)(t _(k) _(j) )′δ_(ij),  (3.4)

where

$\begin{matrix} {{Q_{w}\left( {t_{k_{i}},t_{k_{j}},t_{k_{n}}} \right)} = \begin{Bmatrix} Q_{t_{k_{i}}}^{t_{k_{n}}} & {{{if}\mspace{14mu} t_{k_{i}}} \geq t_{k_{j}}} \\ Q_{t_{k_{j}}}^{t_{k_{n}}} & {{{if}\mspace{14mu} t_{k_{i}}} \geq t_{k_{j}}} \end{Bmatrix}} & (3.5) \\ {{P_{xw}\left( {t_{k_{i}},t_{k_{j}}} \right)} = \begin{Bmatrix} 0 & {{{if}\mspace{14mu} t_{k_{i}}} \leq t_{k_{j}}} \\ {{- {A\left( t_{k_{i}} \right)}}Q_{t_{k_{j}}}^{t_{k_{i}}}{\Phi\left( {t_{k},t_{k_{i}}} \right)}^{\prime}} & {{{if}\mspace{14mu} t_{k_{i} - 1}} < t_{k_{j}} \leq t_{k_{i}}} \\ {{- {A\left( t_{k_{i}} \right)}}Q_{t_{k - 1}}^{t_{k_{i}}}{\Phi\left( {t_{k},t_{k_{i}}} \right)}^{\prime}} & {{{if}\mspace{14mu} t_{k_{j}}} \leq t_{k_{i} - 1}} \end{Bmatrix}} & (3.6) \end{matrix}$ and Q_(c) is defined in eqn. (A.17) in Appendix A.

The solution to the first step is to analytically determine the error covariance matrix.

Theorem 1: The error covariance matrix of the fused track using the fusion rule in eqn. (2.13) is given by

$\begin{matrix} {{P_{f}\left( k \middle| k \right)} = {{\sum\limits_{i = 1}^{n - 1}\;{\sum\limits_{j = 1}^{n - 1}\;{L_{i}M_{ij}L_{j}^{\prime}}}} + {\sum\limits_{i = 1}^{n - 1}\;{L_{i}N_{i}}} + {\sum\limits_{i = 1}^{n - 1}\;{N_{i}L_{i}^{\prime}}} + {M_{n}.}}} & (3.7) \end{matrix}$

For the Proof, see Appendix A.

Additionally, for the actual solution of the problem, the matrices can be defined as

$\begin{matrix} {M = \begin{bmatrix} M_{11} & M_{12} & M_{13} & \cdots & M_{1,{n - 1}} \\ M_{21} & M_{22} & M_{23} & \cdots & M_{2,{n - 1}} \\ \vdots & \; & \ddots & \; & \vdots \\ M_{{n - 1},1} & M_{{n - 1},2} & M_{{n - 1},3} & \cdots & M_{{n - 1},{n - 1}} \end{bmatrix}} & (3.8) \\ {{N = \begin{bmatrix} N_{1} \\ N_{2} \\ \vdots \\ N_{n - 1} \end{bmatrix}},} & (3.9) \end{matrix}$ where matrix M_(ij) and array N_(i) are respectively defined in eqns. (3.1) and (3.3), for 1≦i≦n-1 and 1≦j≦n-1.

Theorem 2: The minimum mean square solution of the asynchronous track fusion problem using the fusion rule in eqn. (2.13) is given by

$\begin{matrix} {{X_{f}\left( t_{k} \middle| t_{k} \right)} = {\sum\limits_{i = 1}^{n}\;{L_{i}{X_{i}\left( t_{k_{i}} \middle| t_{k_{i}} \right)}}}} & (3.10) \\ {{{P_{f}\left( t_{k} \middle| t_{k} \right)} = {{LML}^{~\prime} + {L\; N} + {N^{\prime}L^{\prime}} + M_{n}}},} & (3.11) \end{matrix}$ where

$\begin{matrix} {{\begin{bmatrix} L_{1} & L_{2} & \cdots & L_{n - 1} \end{bmatrix} = {{- N^{\prime}}M^{- 1}\mspace{14mu}{and}}}\;} & (3.12) \\ {L_{n} = {{\Phi\left( {t_{k_{i}},t_{k_{n}}} \right)} - {\sum\limits_{i = 1}^{n - 1}\;{L_{i}{{\Phi\left( {t_{k_{i}},t_{k_{n}}} \right)}.}}}}} & (3.13) \end{matrix}$

For the Proof, see Appendix B.

In the proposed sensor-to-sensor asynchronous track fusion rule, the fusion result at time t_(k-1) is not included in the fusion at time t_(k). Another possible fusion rule can include the previous fusion result as follows:

$\begin{matrix} {{{X_{f}\left( t_{k} \middle| t_{k} \right)} = {\sum\limits_{i = 1}^{n}{L_{i}{X_{i}\left( t_{k_{i}} \middle| t_{k_{i}} \right)}}}},} & (3.14) \end{matrix}$ where for i=0, X ₀(t _(k) ₀ |t _(k) ₀ )=X _(f)(t _(k-1) |t _(k-1)).  (3.15)

In Theorem 1, the sensor observation errors are assumed to be uncorrelated. In cases where the observation noises are coupled across sensors, the only modification needed is in the cross covariance matrix P_(ij). In such cases P_(ij)(t_(ki),t_(kj)) of eqn. (3.4) is replaced with P _(ij)(t _(k) _(i) ,t _(k) _(j) )=A _(i)(t _(k) _(i) )[φ(t _(k) _(i) ,t _(k) _(i) ₋₁)P _(ij)(t _(k) _(i) ₋₁ ,t _(k) _(j) ₋₁)φ(t _(k) _(j) ,t _(k) _(j) ₋₁)+Q _(c)(t _(k) _(i) ,t _(k) _(j) )]A _(j)(t _(k) _(i) )′+K _(i)(t _(k) _(i) )[R _(i)(t _(k) _(i) )δ_(ij) +{hacek over (R)}(t _(k) _(i) ,t _(k) _(j) )]K _(j)(t _(k) _(j) )′  (3.16) where {hacek over (R)}_(ij) (k_(i)) represents the correlation across sensors.

FIG. 1 shows a block diagram of an asynchronous track fusion architecture 100. An instrumentation source set 110 provides an n plurality of multiplexed sources S₁, S₂, S₃, . . . , S_(n). These sources may be identified by the labeled instrumentation S₁ as 111, S₂ as 112, S₃ as 113 and S_(n) as 114. Each source S_(i) provides a corresponding measurement z_(i) and corresponding to a time t_(i) within an interval, where i=1, 2, 3, . . . , n. These measurements provide a measurement set 120 for z_(i) as 121, z₂ as 122, z₃ as 123 and z_(n) as 124.

These measurements are obtained within a time interval between t_(k-1) and t_(k), as shown on a timeline 130. The initial interval time t_(k-1) is marked 131. The first measurement is obtained at time t_(k1) marked 132. The second measurement is obtained at time t_(k2) marked 133. The third measurement is obtained at time t_(k3) marked 134. The n^(th) measurement is obtained at time t_(kn) marked 135. The final interval time t_(k) is marked 136.

A set of local processors 140 provides the continuous time system state X_(i) and the error covariance matrix P_(i) for each measurement, i=1, 2, 3, . . . , n. Note that the local state estimate denoted by X_(i) (t_(ki)|t_(ki)) may be produced by processing the measurement z_(i) (t_(ki)) of sensor i obtained at time t_(ki), and expressed in eqn. (3.10) for the minimum mean square solution. Similarly, the error covariance denoted by P_(i)(t_(ki)|t_(ki)) may be processed by the expression in eqn. (3.11). Each measurement z_(i) and corresponding time t_(ki), is provided with a local processor, such that the first, second, third and n^(th) measurement and time are received in the first, second, third and n^(th) local processors as 141, 142, 143 and 144, respectively. These processors provide a set 150 of the state and covariance values for the first, second, third and n^(th) sets as 151, 152, 153 and 154, respectively. Note that in this embodiment, the local processors implement a Kalman filter to process the local measurements.

Depending on system transmission conditions, the values for the state estimate and error covariance may not arrive sequentially for track fusion. Each value may be delayed by a corresponding lag time λ_(i). For example, a timeline 160 begins at initial interval time t_(k-1) marked 161. The first values may be received at time t_(k1)+λ₁ marked 162; the third values at time t_(k3)+λ₃ marked 163; the second values at t_(k2)+λ₂ marked 164; and the nth value at t_(kn)+λ_(n) marked 165. Note that the second values are received subsequent to the third values due to the much longer second lag time λ₂ as compared to third lag time λ₃. The values are input to an asynchronous track fusion center 170 to yield the fusion values X_(f)(t_(k)|t_(k)) and P_(f)(t_(k)|t_(k)) marked as 180, for merging local values from measurements i=1, 2, 3, . . . , n.

Application -Multitarget-Multisensor Tracking: Multi-target single sensor tracking has developed over many years as sensor capability has expanded. Tracking algorithms have evolved from alpha-beta filters to interactive multiple model filters in more recent years. While the algorithm proposed in this work does not affect single platform multi-target tracking, it does affect target tracking capability when information is being passed back and forth between different platforms in an attempt to synthesize a common air picture, what has been termed by some as a single integrated air picture SIAP. The process incorporating the algorithm described herein represents an enabling technology to synthesize an improved air picture that more closely resembles the operational reality, rather than current technology which involves many over-simplifications. The picture produced by current technology provides a scene that is further in the past than desirable for real-time decision-making. The algorithm presented herein enables a picture to be produced that has fewer latencies and time delays, as well as less uncertainty. The result is a SIAP that is closer to real-time and less “fuzzy” (i.e., more definitive).

Special Case Example—The Bar-Shalon—Compo Fusion Rule: For the case of two asynchronous sensors that arrive at the fusion center, two local estimates can be used to calculate a global estimate. Using Theorem 1 for n=2, the fused track is provided explicitly by X _(f)(t _(k) |t _(k))=L ₁ X ₁(t _(k) ₁ |t _(k) ₁ )+L ₂ X ₂(t _(k) ₂ |t _(k) ₂ )  (3.17) P _(f)(t _(k) |t _(k))=L ₁ M ₁₁ L′ ₁ +L ₁ N ₁ +N′ ₁ L′ ₁ +M ₂,  (3.18) where L ₁ =−N′ ₁ M ₁₁ ⁻¹  (3.19) L ₂=φ(t _(k) ₁ ,t _(k) ₂ )−L ₁φ(t _(k) ₁ ,t _(k) ₂ )  (3.20) M ₁₁ =P ₁₁+φ(t _(k) ₁ ,t _(k) ₂ )[P ₂₂ +Q _(w)(t _(k) ₁ ,t _(k) ₁ ,t _(k) ₂ )]φ(t _(k) ₁ ,t _(k) ₂ )′−P ₁₂(t _(k) ₁ ,t _(k) ₂ )φ(t _(k) ₁ ,t _(k) ₂ )′−φ(t _(k) ₁ ,t _(k) ₂ )P ₂₁(t _(k) ₂ ,t _(k) ₁ ) −P _(xw)(t _(k) ₁ ,t _(k) ₁ )φ(t _(k) ₁ ,t _(k))′−φ(t _(k) ₁ ,t _(k))P _(xw)(t _(k) ₁ ,t _(k) ₁ )′+φ(t _(k) ₁ ,t _(k) ₂ )P(t _(k) ₂ ,t _(k) ₁ )φ(t _(k) ₁ ,t _(k))′+φ(t _(k) ₁ ,t _(k))P(t _(k) ₁ ,t _(k) ₂ )φ(t _(k) ₁ ,t _(k) ₂ )′  (3.21)

$\begin{matrix} {N_{1} = {{{P_{12}\left( {t_{k_{1}},t_{k_{2}}} \right)}{\Phi\left( {t_{k},t_{k_{2}}} \right)}^{\prime}} - {{{\Phi\left( {t_{k_{1}},t_{k_{2}}} \right)}\left\lbrack {{P_{22}\left( {t_{k_{2}},t_{k_{2}}} \right)} - Q_{t_{k_{1}}}^{t_{k_{2}}}} \right\rbrack}{\Phi\left( {t_{k},t_{k_{2}}} \right)}^{\prime}}}} & (3.22) \\ {M_{2} = {{{\Phi\left( {t_{k},t_{k_{2}}} \right)}P_{22}{\Phi^{\prime}\left( {t_{k},t_{k_{2}}} \right)}} + {Q_{t_{k_{2}}}^{t_{k}}.}}} & (3.23) \end{matrix}$

The condition in which the two sensors are synchronous and there are no communication delays between sensor platforms and fusion center, then t_(k1)=t_(k2)=t_(k), and φ(t _(k) _(i) ,t _(k) _(j) )=φ(t _(k) ,t _(k))=1  (3.24)

$\begin{matrix} {Q_{t_{k_{i}}}^{t_{k_{j}}} = {Q_{t_{k}}^{t_{k}} = 0}} & (3.25) \end{matrix}$ P _(xw)(t _(k) ₁ ,t _(k) ₂ )=P _(xw)(t _(k) ,t _(k))=0  (3.26)

Eqns. (3.21)-(3.23) and (3.19)-(3.20) reduce, respectively to M ₁₁ =P ₁₁ +P ₂₂ −P ₁₂ −P ₂₁  (3.27) N ₁ =P ₁₂ −P ₂₂  (3.28) M ₂ =P ₂₂  (3.29) L ₁ =[P ₂₂ −P ₂₁][P ₁₁ +P ₂₂ −P ₁₂ −P ₂₁]⁻¹  (3.30) L ₂ =I−L ₁=(P ₁₁ −P ₁₂)[P ₁₁ +P ₂₂ −P ₁₂ −P ₂₁]⁻¹  ,(3.31) and the fused track eqns. (3.17)-(3.18) becomes X _(f)(k|k)=X ₂(k|k)+(P ₂₂ −P ₂₁)[P ₁₁ +P ₂₂ −P ₁₂ −P ₂₁]⁻¹ [X ₁(k|k)−X ₂(k|k)]  (3.32) P _(f)(k|k)=P ₂₂−(P ₂₂ −P ₂₁)[P ₁₁ +P ₂₂ −P ₁₂ −P ₂₁]⁻¹ [P ₂₂−P₁₂].  (3.33)

Eqns. (3.32)-(3.33) are known as the Bar-Shalom—Campo fusion rule. See, e.g., Y. Bar-Shalom and L. Campo, “The Effect of the Common Process Noise on the Two-Sensor Fused Track covariance”, IEEE Trans. Aerospace & Electronic Systems, vol. AES-22, no. 6, November 1986, pp. 803-805. The case of a plurality of n synchronous sensors can be derived following the same steps.

Concluding Remarks: Various exemplary embodiments provide the solution to a general asynchronous sensor-to-sensor linear distributed estimation problem. On the theoretical level, these embodiments solve three problems in one step: fusing asynchronous tracks without the need for synchronization, dealing with out-of-sequence tracks, as well as accounting for latent tracks. The proposed fusion rule is optimal in the minimum mean square sense.

On the practical level, sensor-to-sensor track fusion can now be performed without the additional computational cost of track synchronization. In addition, using the disclosed solution, the user can include latent and out-of-sequence data in the fusion process without filter re-initialization or back tracking. The rule also applies in the cases when the observation errors are correlated across sensors. In various embodiments, the focus is on sensor-to-sensor track fusion architecture. Other architectures including sensor-to-sensor track fusion with prior and sensor-to-sensor track fusion with feedback can be handled by simple modifications of these embodiments without departing from the spirit of the invention.

Thus, the method and system provided solves an asynchronous sensor-to-sensor track fusion problem using a linear fusion rule. Typical assumptions concerning synchronous operations are removed from the formulation, so the more realistic operational environment is treated in real-time that the sensors data updates are asynchronous. Effectively this means that the algorithm is designed to handle communication delays that exist between sensor platforms and track fusion center, track data can arrive out-of-sequence of the actual time sequence that they were collected. Based on these strenuous real-world operational conditions that occur in the Navy as well as for some type of air platforms, the asynchronous track fusion algorithm described herein deals with all of these contingences.

For the proposed linear fusion rule, the solution is shown to be optimal in the minimum mean square sense. The fusion rule does not necessitate the synchronization of tracks for the purpose of track fusion. The batch processing of incoming tracks provides an elegant solution for the out-of-sequence and latent tracks without special processing of such data. Previous fusion rules, have been demonstrated to be a special case of this rule: the Bar-Shalom—Campo fusion rule that has been commonly applied to two target track fusion problems represents a special case of the asynchronous track fusion algorithm. Furthermore, this algorithm applies not only to track fusion but can also be extended to apply to any other type of state fusion problem.

Appendix A—Proof of Theorem 1: Assume that during the time [t_(k-1), t_(k)], an arbitrary number, n, of asynchronous tracks arrive at the fusion center. Using eqns. (2.18) and (2.19), the error covariance matrix of the fused track is given by

$\begin{matrix} {{P\left( t_{k} \middle| t_{k} \right)} = {\quad{{\sum\limits_{i = 1}^{n - 1}{\sum\limits_{j = 1}^{n - 1}{L_{i}\begin{Bmatrix} {{P_{ij}\left( {t_{k_{i}},t_{k_{j}}} \right)} + {{\Phi\left( {{t_{k}}_{i},t_{k_{n}}} \right)}{P_{nj}\left( {t_{k_{n}},t_{k_{n}}} \right)}{\Phi\left( {t_{k_{j}},t_{k_{n}}} \right)}^{\prime}} +} \\ \begin{matrix} {{{- {P_{i\; n}\left( {t_{k_{i}},t_{k_{n}}} \right)}}{\Phi\left( {t_{k_{j}},t_{k_{n}}} \right)}^{\prime}} - {{\Phi\left( {t_{k_{j}},t_{k_{n}}} \right)}{P\left( {t_{k_{n}},t_{k_{j}}} \right)}} +} \\ \begin{matrix} {{{- {E\left\lbrack {{{\overset{\sim}{X}}_{i}\left( {t_{k_{i}},t_{k_{i}}} \right)}\left( W_{t_{k_{j}}}^{t_{k}} \right)^{\prime}} \right\rbrack}}{\Phi\left( {t,t} \right)}^{\prime}} +} \\ \begin{matrix} {{- {\Phi\left( {t_{k_{i}},t_{k}} \right)}}{E\left\lbrack {{W_{t_{k_{i}}}^{t_{k}}{{\overset{\sim}{X}}_{j}\left( {t_{k_{j}},t_{k_{j}}} \right)}^{\prime}} +} \right.}} \\ \begin{matrix} {{{\Phi\left( {t_{k_{i}},t_{k_{n}}} \right)}{E\left\lbrack {{\overset{\sim}{X}\left( {t_{k_{n}},t_{k_{n}}} \right)}^{\prime}\left( W_{t_{l_{j}}}^{t_{k}} \right)^{\prime}} \right\rbrack}{\Phi\left( {t_{k_{j}},t_{k_{n}}} \right)}^{\prime}} +} \\ \begin{matrix} {{{\Phi\left( {t_{k_{i}},t_{k}} \right)}{E\left\lbrack {\left( W_{t_{k_{i}}}^{t_{k}} \right)^{\prime}{\overset{\sim}{X}\left( {t_{k_{n}},t_{k_{n}}} \right)}^{\prime}} \right\rbrack}{\Phi\left( {t_{k_{j}},t_{k_{n}}} \right)}^{\prime}} +} \\ {{\Phi\left( {t_{k_{i}},t} \right)}{E\left\lbrack {\left( {W_{t_{k_{i}}}^{t_{k}} - W_{t_{k_{n}}}^{t_{k}}} \right)\left( {W_{t_{l_{j}}}^{t_{k}} - W_{t_{k_{n}}}^{t_{k}}} \right)^{\prime}} \right\rbrack}{\Phi\left( {t_{k_{j}},t_{k}} \right)}^{\prime}} \end{matrix} \end{matrix} \end{matrix} \end{matrix} \end{matrix} \end{Bmatrix}L_{i}^{\prime}}}} + {\sum\limits_{i = 1}^{n - 1}{L_{i}\begin{Bmatrix} {{{P\left( {t_{k_{i}},t_{k_{n}}} \right)}{\Phi\left( {t_{k},t_{k_{n}}} \right)}^{\prime}} + {{\Phi\left( {t_{k_{i}},t} \right)}{E\left\lbrack {W_{t_{k_{i}}}^{t_{k}}\left( W_{t_{k_{n}}}^{t_{k}} \right)}^{\prime} \right\rbrack}} +} \\ \begin{matrix} {{{- {\Phi\left( {t_{k_{i}},t_{k_{n}}} \right)}}{P\left( {t,t_{k_{n}}} \right)}{\Phi\left( {t,t_{k_{n}}} \right)}^{\prime}} +} \\ \begin{matrix} {{\left. {\left. {{- {\Phi\left( {t_{k_{i}},t_{k_{n}}} \right)}}E} \right\rbrack W_{t_{k_{i}}}^{t_{k}}{\overset{\sim}{X}\left( t_{k_{n}} \middle| t_{k_{n}} \right)}^{\prime}} \right\rbrack{\Phi\left( {t,t_{k_{n}}} \right)}^{\prime}} +} \\ {{- {\Phi\left( {t_{k_{i}},t_{k}} \right)}}Q_{{\overset{\_}{t}}_{k_{n}}}^{t_{k}}} \end{matrix} \end{matrix} \end{Bmatrix}}} + {\sum\limits_{i = 1}^{n - 1}{\begin{Bmatrix} {{{\Phi\left( {t_{k},t_{k_{n}}} \right)}{P_{ni}\left( {t_{k_{n}},t_{k_{i}}} \right)}} + {{E\left\lbrack {W_{t_{k_{n}}}^{t_{k_{n}}}\left( W_{t_{k_{i}}}^{t_{k}} \right)}^{\prime} \right\rbrack}{\Phi\left( {t_{k_{i}},t_{k}} \right)}^{\prime}} +} \\ \begin{matrix} {{{- {\Phi\left( {t_{k},t_{k_{n}}} \right)}}{P_{nn}\left( {t_{k_{n}},t_{k_{n}}} \right)}{\Phi\left( {t_{k_{i}},t_{k}} \right)}^{\prime}} +} \\ \begin{matrix} {{{- {\Phi\left( {t_{k},t_{k_{n}}} \right)}}{E\left\lbrack {{\overset{\sim}{X}\left( t_{k_{n}} \middle| t_{k_{n}} \right)}\left( W_{t_{k_{i}}}^{t_{k}} \right)^{\prime}} \right\rbrack}{\Phi\left( {t_{k_{i}},t_{k}} \right)}^{\prime}} +} \\ {{- Q_{t_{k_{n}}}^{t_{k}}}{\Phi\left( {t_{k_{i}},t_{k}} \right)}^{\prime}} \end{matrix} \end{matrix} \end{Bmatrix}L_{i}^{\prime}}} + {{\Phi\left( {t_{k},t_{k_{n}}} \right)}{P_{nnn}\left( {t_{k_{n}},t_{k_{n}}} \right)}{\Phi\left( {t_{k},t_{k_{n}}} \right)}^{\prime}} + Q_{t_{k_{n}}}^{t_{k}}}}} & \left( {A{.1}} \right) \end{matrix}$

To simplify the expression for P_(f)(k|k), the expectation terms may be evaluated in eqn. (A.1). Using eqn. (2.6) for t_(ki)≦t_(kj),

$\begin{matrix} \begin{matrix} {W_{t_{k_{i}}}^{t_{k}} = {\int_{t_{k_{i}}}^{t_{k}}{{\Phi\left( {t_{k},\tau} \right)}G{\overset{\_}{W}\ (\tau)}{\mathbb{d}\tau}}}} \\ {= {{\int_{t_{k_{i}}}^{t_{k_{j}}}{{\Phi\left( {t_{k},\tau} \right)}G{\overset{\_}{W}\ (\tau)}{\mathbb{d}\tau}}} + {\int_{t_{k_{j}}}^{t_{k}}{{\Phi\left( {t_{k},\tau} \right)}G{\overset{\_}{W}\ (\tau)}{\mathbb{d}\tau}}}}} \\ {= {{{\Phi\left( {t_{k},t_{k_{j}}} \right)}W_{t_{k_{i}}}^{t_{k_{j}}}} + W_{t_{k_{j}}}^{t_{k}}}} \end{matrix} & \left( {A{.2}} \right) \end{matrix}$ and

$\begin{matrix} {{E\left\{ {W_{t_{k_{i}}}^{t_{k}}\left( W_{t_{k_{j}}}^{t_{k}} \right)}^{\prime} \right\}} = {Q_{t_{k_{j}}}^{t_{k}}.}} & \left( {A{.3}} \right) \end{matrix}$

Similarly, for t_(ki)≧t_(kj), one has

$\begin{matrix} {{E\left\{ {W_{t_{k_{i}}}^{t_{k}}\left( W_{t_{k_{j}}}^{t_{k}} \right)}^{\prime} \right\}} = Q_{t_{k_{j}}}^{t_{k}}} & \left( {A{.4}} \right) \end{matrix}$

The combination of eqns. (A.3) and (A.4) gives

$\begin{matrix} {{E\left\{ {W_{t_{k_{i}}}^{t_{k}}\left( W_{t_{k_{j}}}^{t_{k}} \right)}^{\prime} \right\}} = {{Q_{w}\left( {t_{k_{i}},t_{k_{j}},t_{k}} \right)} = {\begin{Bmatrix} {{Q_{t_{k_{i}}}^{t_{k}}\mspace{11mu}{if}\mspace{14mu} t_{k_{i}}} \geq t_{k_{j}}} \\ {{Q_{t_{k_{j}}}^{t_{k}}\mspace{11mu}{if}\mspace{11mu} t_{k_{i}}} \leq {t_{k}}_{j}} \end{Bmatrix}.}}} & \left( {A{.5}} \right) \end{matrix}$

Using (2.12)

$\begin{matrix} \begin{matrix} {{E\left\{ {{{\overset{\sim}{X}}_{i}\left( t_{k_{i}} \middle| t_{k_{i}} \right)}\left( W_{t_{k_{j}}}^{t_{k}} \right)^{\prime}} \right\}} = {E\left\{ {{A_{i}\left( t_{k_{i}} \right)}\left\lbrack {{{\Phi\left( {t_{k_{i}},t_{k_{i - 1}}} \right)}{{\overset{\sim}{X}}_{i}\left( t_{k_{i} - 1} \middle| t_{k_{i} - 1} \right)}} -} \right.} \right.}} \\ {\left. {\left. W_{t_{k_{i} - 1}}^{t_{k_{i}}} \right\rbrack\left( W_{t_{k_{j}}}^{t_{k}} \right)^{\prime}} \right\} + {E\left\{ {{K_{i}\left( t_{k_{i}} \right)}{V_{i}\left( t_{k_{i}} \right)}\left( W_{t_{k_{j}}}^{t_{k}} \right)^{\prime}} \right\}}} \\ {= {{- {A_{i}\left( t_{k_{i}} \right)}}E{\left\{ \left\lbrack {W_{t_{k_{i} - 1}}^{t_{k_{i}}}\left( W_{t_{k_{j}}}^{t_{k}} \right)}^{\prime} \right\rbrack \right\}.}}} \end{matrix} & \left( {A{.6}} \right) \end{matrix}$

Following the same procedure as in eqn. (A.2), it can be shown that

$\begin{matrix} {{{E\left\{ {{{\overset{\sim}{X}}_{i}\left( {t_{t_{k_{i}}}❘t_{t_{k_{i}}}} \right)}\left( W_{t_{k_{j}}}^{t_{k}} \right)^{\prime}} \right\}} = {P_{xw}\left( {t_{t_{k_{i}}},t_{t_{k_{j}}}} \right)}},} & \left( {A{.7}} \right) \end{matrix}$ where P_(wx) is given in eqn. (3.6). Similarly,

$\begin{matrix} {{{E\left\{ {W_{t_{k_{j}}}^{t_{k}}{X_{j}\left( {t_{k_{j}}❘t_{k_{j}}} \right)}^{\prime}} \right\}} = {P_{wx}\left( {t_{k_{i}},t_{k_{j}}} \right)}},} & \left( {A{.8}} \right) \end{matrix}$ where

$\begin{matrix} {{P_{wx}\left( {t_{k_{i}},t_{k_{j}}} \right)} = {\begin{Bmatrix} 0 & {{{if}\mspace{14mu} t_{k_{j}}} \leq t_{k_{i}}} \\ {{- {\Phi\left( {t_{k},t_{k_{j}}} \right)}}Q_{t_{k_{i}}}^{t_{k_{j}}}{A_{j}\left( t_{k_{j}} \right)}^{\prime}} & {{{if}\mspace{14mu} t_{k_{j} - 1}} \leq t_{k_{i}} \leq t_{k_{j}}} \\ {{- {\Phi\left( {t_{k},t_{k_{j}}} \right)}}Q_{t_{k_{j} - 1}}^{t_{k_{j}}}{A_{j}\left( t_{k_{j}} \right)}^{\prime}} & {{{if}\mspace{14mu} t_{k_{i}}} \geq t_{k_{j} - 1}} \end{Bmatrix} = {{P_{xw}\left( {t_{k_{j}},t_{k_{i}}} \right)}^{\prime}.}}} & \left( {A{.9}} \right) \end{matrix}$

Now, using eqn. (A.5),

$\begin{matrix} {{E\left\lbrack {\left( {W_{t_{k_{i}}}^{t_{k}} - W_{t_{k_{n}}}^{t_{k}}} \right)\left( {W_{t_{k_{j}}}^{t_{k}} - W_{t_{k_{n}}}^{t_{k}}} \right)^{\prime}} \right\rbrack} = {{Q_{w}\left( {t_{k_{i}},t_{k_{k_{j}}},t_{k}} \right)} - {Q_{t_{k_{n}}}^{t_{k}}.}}} & \left( {A{.10}} \right) \end{matrix}$

Using eqn. (2.8)

$\begin{matrix} \begin{matrix} {Q_{t_{k_{j}}}^{t_{k}} = {{\int_{t_{k_{j}}}^{t_{k_{n}}}{{\Phi\left( {t_{k},\tau} \right)}G{\overset{\_}{q}(\tau)}G^{\prime}{\Phi\left( {t_{k},\tau} \right)}^{\prime}\ {\mathbb{d}\tau}}} +}} \\ {\int_{t_{k_{n}}}^{t_{k}}{{\Phi\left( {{t\ }_{k},\tau} \right)}G{\overset{\_}{q}(\tau)}G^{\prime}{\Phi\left( {t_{k},\tau} \right)}^{\prime}{\mathbb{d}\tau}}} \\ {= {{{\Phi\left( {t_{k},t_{k_{n}}} \right)}Q_{t_{k_{j}}}^{t_{k_{n}}}{\Phi\left( {t_{k},t_{k_{n}}} \right)}^{\prime}} + {Q_{t_{k_{n}}}^{t_{k}}.}}} \end{matrix} & \left( {A{.11}} \right) \end{matrix}$

Therefore, for t_(ki)≦t_(kj)

$\begin{matrix} {{E\left\lbrack {\left( {W_{t_{k_{i}}}^{t_{k}} - W_{t_{k_{n}}}^{t_{k}}} \right)\left( {W_{t_{k_{j}}}^{t_{k}} - W_{t_{k_{n}}}^{t_{k}}} \right)^{\prime}} \right\rbrack} = {{\Phi\left( {t_{k},t_{k_{n}}} \right)}Q_{t_{k_{j}}}^{t_{k_{n}}}{{\Phi\left( {t_{k},t_{k_{n}}} \right)}^{\prime}.}}} & \left( {A{.12}} \right) \end{matrix}$

Following the same procedure for t_(ki)>t_(kj), it can be shown that in either case

$\begin{matrix} {{{E\left\lbrack {\left( {W_{t_{k_{i}}}^{t_{k}} - W_{t_{k_{n}}}^{t_{k}}} \right)\left( {W_{t_{k_{j}}}^{t_{k}} - W_{t_{k_{n}}}^{t_{k}}} \right)^{\prime}} \right\rbrack} = {{\Phi\left( {t_{k},t_{k_{n}}} \right)}{Q_{w}\left( {t_{k_{i}},t_{k_{j}},t_{k_{n}}} \right)}{\Phi\left( {t_{k},t_{k_{n}}} \right)}^{\prime}}},} & \left( {A{.13}} \right) \end{matrix}$ where Q_(w) is defined in eqn. (3.5). Applying eqns. (A.7)-(A.8) and (A.12) to (A.1), and using eqns. (3.1)-(3.3), the covariance matrix eqn. (A.1) reduces to

$\begin{matrix} {{P_{f}\left( t_{k} \middle| t_{k} \right)} = {{\sum\limits_{i = l}^{n - 1}{L_{i}M_{ij}L_{j}^{\prime}}} + {\sum\limits_{i = 1}^{n - 1}{L_{i}N_{i}}} + {\sum\limits_{i = 1}^{n - 1}{N_{1}^{\prime}L_{i}^{\prime}}} + {M_{n}.}}} & \left( {A{.14}} \right) \end{matrix}$

In the following, the analytical expression of the cross covariance matrix is derived. The cross covariance matrix P_(ij) between tracks i and j can be defined as

$\begin{matrix} {{P_{ij}\left( {t_{k_{i}},t_{k_{j}}} \right)} = {{E\left\lbrack {{{\overset{\sim}{X}}_{i}\left( {t_{k_{i}}❘t_{k_{i}}} \right)}{{\overset{\sim}{X}}_{j}\left( {t_{k_{j}}❘t_{k_{j}}} \right)}^{\prime}} \right\rbrack}.}} & \left( {A{.15}} \right) \end{matrix}$

The application of (2.12) to evaluate (A.15) yields

$\begin{matrix} \begin{matrix} {{P_{ij}\left( {t_{k_{i}},t_{k_{j}}} \right)} = {{A_{i}\left( t_{k_{i}} \right)}\left\lbrack {{\Phi\left( {t_{k_{i}},t_{k_{i} - 1}} \right)}{P_{ij}\left( {t_{k_{i} - 1},t_{t_{k_{j} - 1}}} \right)}} \right.}} \\ {{\left. {\Phi\left( {t_{k_{j}},t_{t_{k_{j} - 1}}} \right)^{\prime}} \right\rbrack{A_{j}\left( t_{k_{j}} \right)}^{\prime}} +} \\ {{{A_{i}\left( t_{k_{i}} \right)}E\left\{ {W_{t_{k_{i} - 1}}^{t_{k_{i}}}\left( W_{t_{k_{j} - 1}}^{t_{k_{j}}} \right)}^{\prime} \right\}{A_{j}\left( t_{k_{j}} \right)}^{\prime}} +} \\ {K\left( t_{k_{i}} \right)E\left\{ {{V_{i}\left( t_{k_{i}} \right)}{V_{j}\left( t_{k_{j}} \right)}^{\prime}} \right\}{K_{j}\left( t_{k_{j}} \right)}^{\prime}} \\ {= {{A_{i}\left( t_{k_{i}} \right)}\left\lbrack {{{\Phi\left( {t_{k_{i}},t_{k_{i} - 1}} \right)}{P_{ij}\left( {t_{k_{j}},t_{t_{k_{j} - 1}}} \right)}\Phi\left( {t,t} \right)^{\prime}} +} \right.}} \\ {{\left. {Q_{c}\left( {i,j} \right)} \right\rbrack{A_{j}\left( t_{k_{j}} \right)}^{\prime}} +} \\ {{{K_{i}\left( t_{k_{i}} \right)}{R_{i}\left( t_{k_{i}} \right)}{K_{j}\left( t_{k_{j}} \right)}^{\prime}\delta_{t_{k_{i}}t_{k_{j}}}},} \end{matrix} & \left( {A{.16}} \right) \end{matrix}$ where

$\begin{matrix} {{Q_{c}\left( {t_{t_{k_{i}}},t_{t_{j}}} \right)} = {E{\left\{ {W_{t_{k_{i} - 1}}^{t_{k_{i}}}\left( W_{t_{k_{j} - 1}}^{t_{k_{i}}} \right)}^{\prime} \right\}.}}} & \left( {A{.17}} \right) \end{matrix}$

Following the same procedure as in eqn. (A2), it can be shown that

$\begin{matrix} {{Q_{c}\left( {i,j} \right)} = \left\{ \begin{matrix} 0 & {{{if}\mspace{14mu} t_{k_{i}}} \leq {t_{k_{j} - 1}\mspace{11mu}{or}\mspace{14mu} t_{k_{j}}} \leq t_{k_{i} - 1}} \\ {Q_{t_{k_{j} - 1}}^{t_{k_{i}}}{\Phi\left( {t_{k_{j}},t_{k_{i}}} \right)}^{\prime}} & {{{if}\mspace{14mu} t_{k_{i} - 1}} \leq t_{k_{j} - 1} < t_{k_{i}} \leq t_{k_{j}}} \\ {Q_{t_{k_{i} - 1}}^{t_{k_{i}}}{\Phi\left( {t_{k_{j}},t_{k_{i}}} \right)}^{\prime}} & {{{if}\mspace{14mu} t_{k_{j} - 1}} \leq t_{k_{i} - 1} < t_{k_{i}} \leq t_{k_{j}}} \\ {{\Phi\left( {t_{k_{i}},t_{k_{j}}} \right)}Q_{t_{k_{i} - 1}}^{t_{k_{j}}}} & {{{if}\mspace{20mu} t_{k_{j} - 1}} \leq t_{k_{i} - 1} < t_{k_{j}} \leq t_{k_{i}}} \\ {{\Phi\left( {t_{k_{i}},t_{k_{j}}} \right)}Q_{t_{k_{j} - 1}}^{t_{k_{j}}}} & {{{if}\mspace{20mu} t_{k_{i} - 1}} \leq t_{k_{j} - 1} < t_{k_{j}} \leq {t_{k_{i}}.}} \end{matrix} \right.} & \left( {A{.18}} \right) \end{matrix}$

Thus one has fully defined all the terms involved in the computation of the covariance matrix of the fused track. This completes the proof of Theorem 1.

Appendix B—Proof of Theorem 2: Define L=[L ₁ L ₂ . . . L _(n-1)]  (B.1) so using eqns. (A.16), (3.10)-(3.11), and (3.2), one can write P_(f)(t_(k)|t_(k)) as P _(f)(k|k)=LML′+LN+N′L′+M_(n).  (B.2)

Define the performance measure as J=trace[P _(j)(k|k)]=trace[LML′+LN+N′L′+M_(n)].  (B.3)

The weighting matrix L is chosen to minimize the mean square error of the fused track, i.e., to minimize the performance measure J of eqn. (B.3). A necessary condition for optimally is given by

$\begin{matrix} {\frac{\partial J}{\partial L} = {\frac{\partial{{tr}\left\lbrack {P_{f}\left( k \middle| k \right)} \right\rbrack}}{\partial L} = {{{LM} + N^{\prime}} = 0.}}} & \left( {B{.4}} \right) \end{matrix}$

The solution for the weighting matrix L is given by L=−N′M ⁻¹.  (B.5)

While certain features of the embodiments of the invention have been illustrated as described herein, many modifications, substitutions, changes and equivalents will now occur to those skilled in the art. It is, therefore, to be understood that the appended claims are intended to cover all such modifications and changes as fall within the true spirit of the embodiments. 

1. A method for merging data from an n-plurality of multiplexed measurement sources to a decision-maker, comprising: receiving an n-plurality of measurements z_(i) of the data from the respective sources within an acquisition interval [t_(k-1), t_(k)], where i=1, 2, . . . n, wherein each measurement z_(i) corresponds to an associated sampling period t_(ki); using n local processors to process each of the n-plurality of measurements z_(i) to respectively obtain n-pluralities of local state estimates X_(i)(t_(ki)|t_(ki)) and local error covariances P_(i)(t_(ki)|t_(ki)); determining an n-plurality of lag periods λ_(i), each respectively corresponding to a wait duration for obtaining the local state estimates X_(i)(t_(ki)|t_(ki)) and the local error covariances P_(i)(t_(ki)|t_(ki)); offsetting each of the n-plurality of event times t_(ki) by the corresponding each of the n-plurality of the lag periods λ_(i) to provide an n-plurality of offset times t_(ki)+λ_(i); supplying the n-pluralities of the local state estimates X_(i) (t_(ki)|t_(ki)) and the local error covariances P_(i)(t_(ki)|t_(ki)) to a track fusion center in accordance with the n-plurality of offset times; combining the n-plurality of the local state estimates X_(i)(t_(ki)|t_(ki)) as a fusion state estimate X_(f)(t_(k)|t_(k)) and the n-plurality of the local error covariances P_(i)(t_(ki)|t_(ki)) as a fusion error covariance P_(f)(t_(k)|t_(k)) by the track fusion center within the acquisition interval [t_(k-1), t_(k)] for submission to the decision-maker.
 2. The method according to claim 1, wherein the each measurement corresponds to a relation z _(i)(t _(k) _(i) )=H _(i)(t _(k) _(i) )X(t _(k) _(i) )+V _(i)(t _(k) _(i) ), where H_(i)(t_(ki)) is a measurement matrix of a source i in the n-plurality of sources, X(t_(ki)) is a system state at sampling period t_(ki), and V_(i)(t_(ki)) is measurement noise of the source i.
 3. The method according to claim 1, wherein summing the n-plurality of local state estimates X_(i)(t_(ki)|t_(ki)) further includes: calculating an n-plurality of weighting factors L_(i), each respectively corresponding to the local state estimates X_(i)(t_(ki)|t_(ki)); multiplying the each of the n-plurality of local state estimates X_(i)(t_(ki)|t_(ki)) with a respective each of the n-plurality of weighting factors L_(i) to produce an n-plurality of weighted state estimates L_(i)X_(i)(t_(ki)|t_(ki)); and summing the n-plurality of weighted state estimates L_(i)X_(i)(t_(ki)|t_(ki)) as the fusion state estimate X_(f)(t_(k)|t_(k)).
 4. The method according to claim 1, wherein summing the n-plurality of the local error covariances P_(i)(t_(ki)|t_(ki)) as the fusion error covariance P_(f)(t_(k)|t_(k)) includes minimizing a mean square error term.
 5. The method according to claim 4, wherein the fusion error covariance P_(f)(t_(k)|t_(k)) is determined by a relation ${{P_{f}\left( k \middle| k \right)} = {{\sum\limits_{i = 1}^{n - 1}{\sum\limits_{j = 1}^{n - 1}{L_{i}M_{ij}L_{j}^{\prime}}}} + {\sum\limits_{i = 1}^{n - 1}{L_{i}N_{i}}} + {\sum\limits_{i = 1}^{n - 1}{N_{i}L_{i}^{\prime}}} + M_{n}}},$ where prime symbol ′ represents transpose, M is an error difference matrix and N is an error difference array, such that each element therein is defined by M _(ij) =P _(ij)(t _(k) _(i) ,t _(k) _(j) )+φ(t _(k) _(i) ,t _(k) _(n)) [P _(nn)(t _(k) _(n) ,t _(k) _(n) )+Q _(w)(t _(k) _(i) ,t _(k) _(j) ,t _(k) _(n) )]φ(t _(k) _(j) ,t _(k) _(n) )′−P _(in)(t _(k) _(i) ,t _(k) _(n) )φ(t _(k) _(j) ,t _(k) _(n) )′−φ(t _(k) _(i) ,t _(k) _(n) )P _(nj)(t _(k) _(n) ,t _(k) _(j) ) −P _(xw)(t _(k) _(i) ,t _(k) _(j) )φ(t _(k) _(j) ,t _(k))′−φ(t _(k) _(i) ,t _(k))P _(xw)(t _(k) _(j) ,t _(k))′+φ(t _(k) _(i) ,t _(k) _(n) )P _(xw)(t _(k) _(n) ,t _(k) _(j) )φ(t _(k) _(j) ,t _(k))′+φ(t _(k) _(i) ,t _(k))P _(xw)(t _(k) _(n) ,t _(k) _(i) )′φ(t _(k) _(j) ,t _(k) _(n) )′ P_(ij)(t _(k) _(i) ,t _(k) _(j) )=A(t)[φ(t _(k) _(i) ,t _(k) _(i) ₋₁)P(t _(k) _(i) ₋₁ ,t _(k) _(j) ₋₁)φ(t _(k) _(j) ,t _(k) _(j) ₋₁)′+Q _(c)(t _(k) _(i) ,t _(k) _(j) )]A _(j)(t _(k) _(j) )′+K _(i)(t _(k) _(i) )R _(i)(t _(k) _(i) )K _(j)(t _(k) _(i) )′δ_(if) for 1≦i≦n−1 and 1≦j≦n−1, wherein ${Q_{w}\left( {t_{k_{i}},t_{k_{j}},t_{k_{n}}} \right)} = \begin{Bmatrix} Q_{t_{k_{i}}}^{t_{k_{n}}} & {{{if}\mspace{14mu} t_{k_{i}}} \geq t_{k_{j}}} \\ Q_{t_{k_{j}}}^{t_{k_{n}}} & {{{if}\mspace{14mu} t_{k_{i}}} \leq t_{k_{j}}} \end{Bmatrix}$ ${P_{xw}\left( {t_{k_{i}},t_{k_{j}}} \right)} = {\begin{Bmatrix} 0 & {{{if}\mspace{14mu} t_{k_{i}}} \leq t_{k_{j}}} \\ {{- {A\left( t_{k_{i}} \right)}}Q_{t_{k_{j}}}^{t_{k_{i}}}{\Phi\left( {t_{k},t_{k_{i}}} \right)}^{\prime}} & {{{if}\mspace{14mu} t_{k_{i} - 1}} < t_{k_{j}} \leq t_{k_{i}}} \\ {{- {A\left( t_{k_{i}} \right)}}Q_{t_{k_{i} - 1}}^{t_{k_{i}}}{\Phi\left( {t_{k},t_{k_{i}}} \right)}^{\prime}} & {{{if}\mspace{14mu} t_{k_{j}}} \leq t_{k_{i} - 1}} \end{Bmatrix}{\mspace{11mu}\;}{and}}$ Q_(c)(t_(t_(k_(i))), t_(t_(j))) = E{W_(t_(k_(i − 1)))^(t_(k_(i)))(W_(t_(k_(j) − 1))^(t_(k_(j))))^(′)}, and such that φ(t_(k) _(i) ,t_(k) _(i) ₋₁)=e ^(A(t) ^(ki) ^(−t) ^(ki−1) ⁾ is a transition matrix for a filter difference matrix A over a time interval (t_(ki), t_(ki-1)), W is state process noise, E is an expected value, Q is a process noise covariance matrix, wherein the filter difference matrix A is expressed as A_(i)(t_(k) _(i) )=[I−K_(i)(t_(k))H_(i)(t_(k) _(i) )], such that I is an identity matrix, K_(i)(t_(ki)) is a process filter gain and H_(i)(t_(ki)) is a measurement matrix of a source i of the n-plurality of sources.
 6. The method according to claim 1, wherein summing the n-plurality of the local state estimates X_(i)(t_(k) _(i) |t_(ki)) as the fusion state estimate X_(f)(t_(k)|t_(k)) and the n-plurality of the local error covariances P_(i)(t_(ki)|t_(ki)) as the fusion error covariance P_(f)(t_(k)|t_(k)) are minimum mean square solutions of respective relations ${X_{f}\left( t_{k} \middle| t_{k} \right)} = {\sum\limits_{i = 1}^{n}{L_{i}{X_{i}\left( t_{k_{i}} \middle| t_{k_{i}} \right)}}}$ P_(f)(t_(k)|t_(k)) = L M L^(′) + L N + N^(′)L^(′) + M_(n), wherein $\begin{bmatrix} L_{1} & L_{2} & \ldots & L_{n - 1} \end{bmatrix} = {{- N^{\prime}}M^{- 1}\mspace{20mu}{and}}$ ${L_{n} = {{\Phi\left( {t_{k_{i}},t_{k_{n}}} \right)} - {\sum\limits_{i = 1}^{n - 1}{L_{i}{\Phi\left( {t_{k_{i}},t_{k_{n}}} \right)}}}}},$ L is an array of weighting factors, M is an error difference matrix, N is an error difference array, prime symbol ′ is a transpose and φ is a transition matrix.
 7. The method according to claim 1, wherein the measurements z_(i) correspond to translational parameters.
 8. The method according to claim 7, wherein the translational parameters comprise at least one of range, elevation and azimuth.
 9. The method according to claim 8, wherein the fusion state estimate X_(f) comprises at least one of a position and a velocity with respect to each of the n-plurality of sources i. 